#include "prot-filer/backbones_utils.h"

namespace backbones_utils
{
using namespace mili;

void semilla(BasicProtein& atm, float* R)
{
    /* Los valores semilla */
    atm[0].x = 0.000;
    atm[0].y = 0.000;
    atm[0].z = 0.000;
    atm[1].x = 1.460;
    atm[1].y = 0.000;
    atm[1].z = 0.000;
    atm[2].x = 2.011;
    atm[2].y = 1.422;
    atm[2].z = 0.000;
    atm[0].vdw = N  ;
    atm[1].vdw = CA ;
    atm[2].vdw = C;


    R[0]  = 0.361594;
    R[1]  = -0.932336;
    R[2]  = 0.000000;
    R[3]  = 2.011431;
    R[4]  = 0.932336;
    R[5]  = 0.361594;
    R[6]  = 0.000000;
    R[7]  = 1.421812;
    R[8]  = 0.000000;
    R[9]  = 0.000000;
    R[10] = 1.000000;
    R[11] = 0.000000;
    R[12] = 0.000000;
    R[13] = 0.000000;
    R[14] = 0.000000;
    R[15] = 1.000000;

}

void angles2coord(BasicProtein& protein, const AnglesMapping& angles_mapping, const AnglesData& angles_data)
{
    // This is also present in petu.cpp. Should be hardcoded?
    // Maximun gyration radius and maximun CA-CA distance.
    // Both equations constructed from database analisys.
    // rgmax not necessary.
    //float dmax2 = (8.0 * mili::cubic_root(float(angles_data.nres)) + 25.0);
    float R_inicial[16];
    semilla(protein, R_inicial);
    DummyFilter filter;
    for (unsigned int i = 2; i <= angles_data.nres; ++i)
    {
        const AngleValue fi = angles_mapping.get_fi_value(angles_data.angles[i - 2].fi);
        const AngleValue si = angles_mapping.get_si_value(angles_data.angles[i - 2].si);
        poneres(R_inicial,
                cos(deg2rad(si)),
                sin(deg2rad(si)),
                cos(deg2rad(fi)),
                sin(deg2rad(fi)),
                protein,
                i,
                filter);
    }
}

void int2car(float* pT, float bon, float cosang, float sinang,
             float costor, float sintor, BasicProtein& patm, int at, AtomType tipo)
{
    float A[16], rot[11];

    rot[0] = -cosang;
    rot[1] = -sinang;
    rot[2] =  costor;
    rot[3] = -sintor;
    rot[4] =  bon * rot[0];
    rot[5] = -rot[1] * rot[2];
    rot[6] =  rot[0] * rot[2];
    rot[7] =  bon * rot[5];
    rot[8] =  rot[1] * rot[3];
    rot[9] = -rot[0] * rot[3];
    rot[10] = bon * rot[8];

    A[0]  = pT[0] * rot[0] + pT[1] * rot[5] + pT[2] * rot[8];
    A[4]  = pT[4] * rot[0] + pT[5] * rot[5] + pT[6] * rot[8];
    A[8]  = pT[8] * rot[0] + pT[9] * rot[5] + pT[10] * rot[8];
    A[12] = pT[12] * rot[0] + pT[13] * rot[5] + pT[14] * rot[8];

    A[1]  = pT[0] * rot[1] + pT[1] * rot[6] + pT[2] * rot[9];
    A[5]  = pT[4] * rot[1] + pT[5] * rot[6] + pT[6] * rot[9];
    A[9]  = pT[8] * rot[1] + pT[9] * rot[6] + pT[10] * rot[9];
    A[13] = pT[12] * rot[1] + pT[13] * rot[6] + pT[14] * rot[9];

    A[2]  = pT[1] * rot[3] + pT[2] * rot[2];
    A[6]  = pT[5] * rot[3] + pT[6] * rot[2];
    A[10] = pT[9] * rot[3] + pT[10] * rot[2];
    A[14] = pT[13] * rot[3] + pT[14] * rot[2];

    A[3]  = pT[0] * rot[4] + pT[1] * rot[7] + pT[2] * rot[10] + pT[3];
    A[7]  = pT[4] * rot[4] + pT[5] * rot[7] + pT[6] * rot[10] + pT[7];
    A[11] = pT[8] * rot[4] + pT[9] * rot[7] + pT[10] * rot[10] + pT[11];
    A[15] = pT[12] * rot[4] + pT[13] * rot[7] + pT[14] * rot[10] + pT[15];

    patm[at].x   = A[3];
    patm[at].y   = A[7];
    patm[at].z   = A[11];
    patm[at].vdw = tipo;

    copymat(pT, A);
}
}
